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Abstract 

One goal of this paper is to discuss the classical definition of granular temperature as 
an extension of its thermodynamic equivalent and a useful concept which provides an 
important characterization of fluidized granular matter. Following a review of some 
basic concepts and techniques (with emphasis on fundamental issues) we present 
new results for a system that can exhibit strong violations of equipartition, yet 
is amenable to description by classical granular hydrodynamics, namely a binary 
granular gas mixture. A second goal of this article is to present a result that pertains 
to dense granular and molecular solids alike, namely the existence of a correction 
to the elastic energy which is related to the heat flux in the equations of continuum 
mechanics. The latter is of the same (second) order in the strain as the elastic energy. 
Although recent definitions of temperatures for granular matter, glasses and other 
disordered many-body systems are not within the scope of this article we do make 
several general comments on this subject in the closing section. 
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1 Introduction 



The axiomatic formulation of thermodynamics provides a precise definition of 
temperature. Nearly any textbook on thermodynamics shows this definition 
to be compatible with common experience. Equilibrium statistical mechanics 
endows temperature with a microscopic meaning but invokes a few assump- 
tions (e.g., ergodicity). The physical relevance and usefulness of both of these 
theories is beyond debate. 

As there are practically no equilibrium systems in nature there is clearly a 
need to describe non-equilibrium systems. The most successful theories that 
treat such systems are those devoted to near-equilibrium states, namely states 
in which the macroscopic fields vary sufficiently slowly in space and time so 
that (near) local equilibrium distributions can establish themselves everywhere 
while the state of the system hardly changes. That this should occur is not 
entirely obvious. It helps to note that e.g., in gases local equilibration occurs 
within a few collisions per particle, i.e. on time scales that are much smaller 
than macroscopic times. In other words, these approaches are effective when 
there is significant scale separation between the spatial and temporal micro- 
and-macroscales. In the latter case one can justify hydrodynamics, elasticity 
and other macroscopic theories, and the transport coefficients are given by 
Green-Kubo expressions derived using linear (and non-linear) response theo- 
ries (1) or projection operator approaches (2; 3), a similar statement holding 
for fluctuation dissipation 'theorems' e.g., for Brownian motion (4). 

The macroscopic fields that represent the states of systems that possess scale 
separation are usually taken to be the densities of the conserved entities, 
namely the (mass or number) density, the momentum density and the en- 
ergy density. These fields are "slow" since any change in their local values 
necessitates motion to (non-locally) transport an amount of the entity repre- 
sented by the field. It is usually assumed that on macroscopic time scales the 
system 'remembers' only the slow fields, the other degrees of freedom rapidly 
accommodating themselves to those fields and becoming "noise" whose statis- 
tics is determined by the slow fields. The description of far-from-equilibrium 
systems, such as turbulent fiuids and glasses, is not as intuitively founded as 
that of near-equilibrium systems; in particular, since many of these systems 
lack scale separation it is a-priori unclear which fields are sufficient to pro- 
vide a closed description for them (nor is the nature of a relevant statistical 
ensemble for the unresolved degrees of freedom obvious). 

The BBGKY hierarchy (see e.g., (5)) is an exact formulation for many-body 
systems that does not require them to be reversible (6) on the microscopic 
(particle scale) level or characterized by 'slow' or hydrodynamic fields. How- 
ever, without approximations that truncate it to a manageable size it is nearly 
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useless. A truncation which is vahd at low densities (gases) begets the Boltz- 
mann equation (5; 7). Like the BBGKY hierarchy, the Boltzmann equation 
is well defined without the introduction of the concept of temperature or 
other macroscopic fields and it resolves near-microscopic scales. However it is 
a nontrivial integrodifferential equation that can be directly solved in very few 
cases. One of the important uses of the Boltzmann equation is the derivation 
of hydrodynamic constitutive relations. The form of the hydrodynamic equa- 
tions (i.e., the equations of continuum mechanics) is very general and holds 
for single realizations (see Appendix A). They also follow from the Boltzmann 
equation (whose unknown is a distribution function). A gradient (Chapman- 
Enskog) expansion of the solution of the Boltzmann equation begets explicit 
expressions for the hydrodynamic fiuxes in terms of gradients of the hydro- 
dynamic fields (constitutive relations). This expansion is formally valid when 
there is scale separation (see more below) although the resulting hydrodynamic 
equations have been successfully employed when scale separation is weak or 
non-existent, such as in the case of weak shocks or granular matter (8). 

Granular systems are different from molecular systems not only in the re- 
spective sizes of their constituents but (more importantly) the fact that the 
grains interact in a dissipative way (e.g., collide inelastically in the granular 
gas). Therefore all granular states are of non-equilibrium nature. An impor- 
tant consequence of this fact (9) is that granular systems do not possess strong 
scale separation except in the near-elastic case. In addition, typical granular 
systems comprise much fewer than (say) an Avogadro number of grains, there- 
fore they are very small by molecular standards. Therefore granular systems 
are basically mesoscopic (9): they are expected to exhibit strong fluctuations 
within any reasonably deflned ensemble and realizations may not be represen- 
tative of these ensembles, see e.g., (10; 11), even in the fluidized phase. 

Imagine a granular gas which so close to being elastic that the percentage of 
the energy lost in a collision is minute. Since local equilibration is established in 
a matter of a few collisions per particle the main effect of the weak inelasticity 
is to cause the decay of the kinetic energy on long time scales (of course, if the 
system is forced this 'lost' energy can be replenished). Therefore near-elastic 
granular gases are expected to behave in a way that is close to that of elastic 
gases, as detailed calculations show (8). In other words, the elastic limit is not 
singular (for finite times) and one expects the temperature to be a relevant 
field with properties that are similar to those encountered in elastic systems. 
It is a-priori unclear whether a similar statement holds for strongly inelastic 
granular gases. It is interesting to note that even in the realm of near-elastic 
granular gases there is no known equivalent of Boltzmann's H-theorem. 

Granular solids are jammed in static positions with no possibility to ex- 
plore phase space. They are athermal. Furthermore, their states are typically 
metastable (the ground state of a sand-pile is one in which all grains reside 
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on the ground). Although the concept of temperature seems to be a-priori 
irrelevant to these systems, there have been interesting suggestions to define 
effective (configurational) temperatures for them as properties of appropriate 
ensembles, e.g., (12; 13). Slow or quasistatic flows of dense granular matter 
are clearly athermal as well. Here possible dynamic applications of the above 
mentioned effective temperatures, such as fluctuation-dissipation (FD) rela- 
tions have been proposed, and tested, see e.g., (14; 15), indicating a possible 
dynamical relevance of these temperatures. Unlike the FD relations for gran- 
ular gases (16; 17; 18) (which are not trivial but can be derived by means 
borrowed from near-equilibrium statistical mechanics), the FD relations for 
dense and quasistatic granular flows (or glasses) have been conjectured (see, 
though (19)). Still the above mentioned as well as other flndings are interesting 
and suggestive of a possible underlying theoretical framework. Such a frame- 
work is probably not going to enjoy the generality of the near-equilibrium 
formulations and may necessitate an increase in the number of variables that 
characterize the considered systems (as is often the case in non-equilibrium 
situations). A proposal for deflning non-equilibrium temperatures which has 
attracted considerable attention is due to Tsallis (20). Whether this deflnition 
is relevant to granular matter is unclear as yet, see though (21). 

In old theories of granular matter (gases included) the granular temperature 
was not accounted for (as it does not correspond to a conserved entity) (22). 
Even in some modern studies the granular temperature is (used but) elim- 
inated from the equations of motion [e.g., for length scales for which heat 
conduction is subdominant to heat production and dissipation (23)]; in these 
cases the temperature is claimed to be enslaved to the other flelds (but it still 
plays a useful role as a state characterizing fleld). Systematic derivations of 
hydrodynamic equations for granular gases (24; 25) invariably account for the 
granular temperature as a relevant fleld. These equations have been validated 
against simulations and experiments (8), even when scale separation is weak. 

The present paper has three major goals. The flrst is to show that the temper- 
ature fleld is needed for the description of granular g clSGS clS db characterizing 
fleld, if not a quasi-thermodynamic entity. This is done in Section 2, where the 
concept of granular temperature is explained in the context of the CE expan- 
sion. The goal of Section 3 is to provide an example of a kinetic theory (for 
binary granular gases) in which the temperature fleld is one of the hydrody- 
namic flelds and which does produce sensible results even when equipartition 
is strongly violated (here a novel method for carrying out the CE expansion 
is presented and applied). The third goal, taken up in Section 4, is to show 
the existence of a correction to the classical formula for the elastic energy of 
granular and molecular matter alike; this correction is shown to be related to 
the heat flux in the equations of continuum mechanics. Section 5 concludes 
this paper and briefly discusses some issues concerning non-thermodynamic 
temperatures. Some technical details are relegated to Appendices. 
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2 Granular Gases and the Chapman-Enskog expansion 



As mentioned in the Introduction, the equations of continuum mechanics are 
very general and in particular, they hold for practically all single realizations 
of many body systems (see Appendix A). These equations comprise the equa- 
tion of continuity (for the density, p), Eq. (A. 2), the momentum density, p, 
equation, Eq. (A. 10), and the equation for the energy density, e, Eq. (A. 15), 
in which (in the dissipative/inelastic case) an energy sink term, pF, arises. 
An equation for the internal specific energy easily follows from the above, see 
Eq. (A. 20). Note that the derivation in Appendix A produces also explicit 
expressions for the relevant fields and fluxes in terms of the microscopic de- 
grees of freedom. In particular, note that the expression for the stress tensor 
is not precisely the same as the commonly used Voigt or Born or Kirkwood 
expression; the latter are limits of Eq. (A. 9) when the coarse-graining volume 
(or width of the coarse-graining function) is very large and do not satisfy 
the equations of continuum mechanics unless there is scale separation or one 
considers homogeneous deformations of lattice configurations. 

In the realm of molecular or granular gases, when the internal potential energy 
is negligible (for dilute systems) or when a model of hard spheres is employed 
one can replace the specific internal energy by the (granular) temperature, T. 
In this case Eq. (A. 20) can be rewritten as: 



p(r,t)^^ = ^^lj^ap^{r,t) - div Q(r, t) - p(r, t)r(r, t) (1) 



with notation that is defined in Appendix A. Following Eq. (A. 19) the (gran- 
ular) temperature field is given by the following intuitive expression: 



Note that when scale separation does not exist, the temperature, much like 
other fields, may be resolution dependent (i.e. it depends on the width of the 
coarse graining function 0)(26). 

The hydrodynamic description of a fluid, granular or otherwise, is complete 
when constitutive relations express cr, Q and F in terms of the hydrodynamic 
fields. It is a-priori unclear whether such a closure exists for granular fluids. 
Extended Chapman-Enskog (CE) expansions provide explicit constitutive re- 
lations for granular g further detailed below. 

The Boltzmann equation is a non-trivial integrodifferential equation for the 
single-particle velocity distribution, /(v,r,t), whose integral over the (parti- 
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cle) velocities, v, is the number density, n(r, t). The ratio /(v, r, t)/?i(r, t) is 
the probabihty density for a particle to have velocity v when it is found (or 
more accurately its center of mass is found) at point r at time t. In the rela- 
tively simple case of a dilute monodisperse collection of smooth hard spheres 
of unit mass and diameter d whose collisions are characterized by a fixed co- 
efficient of normal restitution, a, the Boltzmann equation reads (24; 25; 27): 

= d' J rfV2rfk(k ■ V12) (^/(V;)/(V^) - /(Vi)/(V2)) (3) 
kvi2>0 

where V is a spatial gradient. The unit vector k points from the center of 
particle T to the center of particle '2'. The relation between the precollisional 
(primed) and postcoUisional velocities is given by: Vi = v- — ^-|^(k ■ V;j)k 
where v-j = v( — Vj. The dependence of / on the spatial coordinates and 
time is not explicitly spelled out in Eq. (3). Notice that in addition to the 
exphcit dependence of Eq. (3) on a, it also implicitly depends on a through the 
relation between the postcoUisional and precollisional velocities. The condition 
k ■ V12 > represents the fact that only particles whose relative velocity is 
such that they approach each other, can collide. 

In the derivation of the Boltzmann equation (5; 7) it is assumed that the veloc- 
ities and positions of colliding particles are not correlated. This assumption 
of "molecular chaos" is justified for dilute gases but not for dense gases or 
strongly inelastic granular gases (9; 28; 29; 30). This fact is related to the re- 
duction in the normal component of the relative velocity of colliding particles. 
In relatively dense gases particles heading towards each other may not collide 
as they are 'stopped' by collisions with other particles. This phenomenon is 
modeled in the Enskog-corrected Boltzmann equation (7). 

The Chapman- Enskog expansion is a method to perturbatively solve the Boltz- 
mann equation. Scale separation is formally necessary for its implementation. 
When all gradients vanish, a local Maxwellian distribution, corresponding to 
the values of the hydrodynamic fields at each point, solves the Boltzmann 
equation in the elastic case. A gradient expansion (formally, an expansion in 
the Knudsen number) around this solution automatically yields the result that 
the distribution function depends on space and time only through its depen- 
dence on the fields (which is quasi-local since f{r,t) depends on the values 
of the fields and their spatial derivatives at (r, t)). The often made statement 
that within the CE expansion / is assumed to depend on space and time 
only through its dependence on the hydrodynamic fields is not precise: this 
assumption is made only for the zeroth order of the expansion. 



6 



In the case of granular gases, the (kinetic) energy density is not strictly con- 
served, and this led to questions whether it should be included in the set 
of hydrodynamic fields (8; 23; 31; 32). Following the discussion in the intro- 
duction, when the degree of inelasticity is sufficiently small, it is justified to 
include the granular temperature in the set of hydrodynamic fields. An ex- 
tension of the Chapman Enskog expansion based on this approach has been 
proposed and implemented in (24). There, the small parameters used to ex- 
pand the Boltzmann equation were the Knudsen number, K = where i is 
the mean free path given by i = and L is a macroscopic length scale 

(i.e., the length scale which is resolved by hydrodynamics, not necessarily the 
system size) and the degree of inelasticity, defined by e = 1 — a^. In another 
approach (25) to the problem only the Knudsen number is the small param- 
eter and the temperature is assumed to comprise a relevant hydrodynamic 
field for all values of the degree of inelasticity. This method does not seem to 
lead to contradictions for any value of a; as a matter of fact it has recently 
been used to obtain accurate values for the linear transport coefficients for all 
physically allowed values of a (33). Whether the results are indeed relevant 
to strongly inelastic granular gases, where scale separation does not hold and 
precoUisional correlations may be prominent, remains to be seen. 



3 Binary granular gases: extreme violations of equipartition 



Consider a mixture of smooth hard spheres, composed of species A and B, of 
masses and m^, and diameters a a and a^, respectively. The coefficient of 
normal restitution (assumed to be fixed) for a collision of a particle of species 
a G {A, B} with a particle of species /3 G {A, B} is denoted by eap (hence, 
60/3 G {caa, ^bb, ^ab})- The transformation of velocities due to a collision of 
a sphere of species a with a sphere of species P is given by: 



vi = v;-(l + e,^)M"'3(v;2-A;)A; (4) 
V2 = + (1 + e^p) M-^ (v;2 ■ k) k (5) 

where {v']^,V2} denote the precoUisional velocities of the spheres (the index 
'1' refers here to species a), and {vi, V2} are the corresponding postcollisional 
velocities; A; is a unit vector pointing from the center of sphere a to that of 
sphere /5, and M"^ = — ^ — V12 = vi — V2, a similar definition holding for 
the primed (precoUisional) velocities. Obvious kinematic constraints require 
that v'^2 ■ k > 0. Define the degrees of inelasticity as: eai3 = 1 — e^^. 

The kinetic description of a binary granular gas mixture involves two Boltz- 
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mann equations, one for each species (34; 35): 



dt 



+ ^1 -Vfa = B{faja,eaa) + B{faJ(3,eaf3) a ^ (3, (6) 



where /„ (v) is the distribution function for particles of species a, and 



vi2k>0 



/«(v;)//3K) 



- /a(vi)//3(v2) 



V12 ■ k (iv2 dk, 



(7) 



where cXo/? = 



Here vi and v'^ pertain to species a, and V2 and 
pertain to (3. Notice that B{fa, f/3, eap) depends on the coefficients of normal 
restitution both explicitly, as shown in Eq. (7), and implicitly through the 
collision law. 

The hydrodynamic fields in the present case comprise the two number den- 
sities. Ha and ub (or the mass densities pa = nA^A and ps = nBrris), the 
mixture's velocity field, V, and the temperature field, T, defined (differently 
from the monodisperse case) as twice the mean fluctuating kinetic energy of a 
particle (these are the only slow flelds in the elastic limit), as explicitly deflned 
below. The continuum equations of motion follow from the pertinent Boltz- 
mann equation directly (using the standard procedure of computing velocity 
moments of the equation). Their validity is general since they are based on the 
underlying conservation laws. The equation of motion for the number density, 
Ua, with a G {A,B}, is: 



Dt 



— div J„ — n„ div V 



where is the material derivative and Jq, = Ua (Va — V) is the particle flux 
density of species a. As Vq,, the velocity fleld of species a (or the flux Jq,) 
is not a hydrodynamic fleld, it must be given by an appropriate constitutive 
relation. The velocity fleld obeys (as expected): = — + Pg? where 

the summation convention is used, and 0"^ is the stress tensor. The granular 
temperature fleld obeys: n-^ = T div J — div Q — 2aij^ — F, where 
J = + is the total particle flux, and Q is the heat flux. 

The form of the constitutive relations can be easily determined from tenso- 
rial (and symmetry) considerations, the result being (to linear order in the 
gradients, or Navier-Stokes order): aij = p6ij — 2p Dij — 6ij rjj^ div V, where 

~ I ( §^ " I divV j is the traceless rate of strain tensor, p is 
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the shear viscosity, ?7b is the bulk viscosity (which vanishes in the dilute limit) 
and p is the pressure. The diffusion flux can be rewritten as 



n a 



AB 



<Vlnn-<Vlnc 



(9) 



where k^, and are non-trivial functions of the parameters 

S = {{^a/?}, Ma = ^ , ^ , c}, 

where the concentration field, c, is defined as c = — = — — . Similarly, the 
heat flux can be rewritten as 

1 T^/^ / ^ \ 
Q = — - A^VhiT- A"Vlnn- A^Vlnc (10) 



where mo = + m^. In the dilute limit the equation of state is the same 
as that for an ideal gas: p = Much like in the monodisperse case the 
number density for species a is given by: Ua = J /a(v) dv, the corresponding 
mass density being = maUa] the overall number density is n = + n^, 
and the overall mass density is p = + Pb- The velocity field of species 
a is given by: Vq, = — / /a(v) v rfv. As mentioned, is not a hydrody- 
namic field, and needs to be expressed as a functional of the hydrodynamic 
fields. The mixture's velocity field is: V = ^ {pa^a + Pb^b)- The granular 
temperature of species a is defined by: Ta = ^ faiy) f^a (v — V)^ c/v 
(not a hydrodynamic field). The velocity fluctuations of each of the species 
are measured with respect to the (hydrodynamic) mixture's velocity field, not 
the species' velocity fields. The mixture's granular temperature is defined as: 
T = ^ {uaTa + ubTb)- The kinetic expression for the stress tensor (7) is: 
(Tjj = rriA J /a (v) UiUj dv + niB J fs (v) UiUj c/v, where u = v — V is the 
peculiar (fluctuating) velocity of a particle (irrespective of the species). Sim- 
ilarly, the heat flux is composed of two contributions, for obvious reasons: 
Q = / /a (v) rriAU^u dv + J /b (v) itibu^u dv. Our definition of the heat flux 
may differ by a factor of 2 from some other definitions. The sink term is given 
by r = + + r^^, where: 

2 



and 

TaB^Sab"^^^^^^^ J J /A(Vl)/,,(v2)|Vi2|'rfVi rfV2. (12) 
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The HCS distributions for binary granular gas mixtures (much hke in the 
monodisperse case) are scahng solutions, in which the distribution functions 
are rendered time independent when the velocities are scaled by the respective 
(decaying) thermal speeds of the species comprising the mixture. It serves as 
a zeroth order in the extended Chapman- Enskog expansion of (34). The Ihs 
of the Boltzmann equations for the HCS can be rewritten in this case as: 

Qf HCS Q f HCS r,rp pQJ-HCS 

= ^ = — - grp 0.10.0. t hc two Boltzmauu equations become: 



n dT 



B /r^ c'e- + B {fT\ fr^a,) (is) 



These equations need to be solved subject to the following 'constraints' which 
fix the values of the hydrodynamic fields: / fj^^du = ua, J fs'^^du = ub 
and / fA^^mAu'^du+ J fg'^^mBU^dvL = nT. The common way of solving Eqs. 
(13) is to represent the functions, f^^^, by Sonine polynomial series times 
Maxwellians in the respective non-dimensionalized (peculiar) velocities. How- 
ever, it can be shown (33) that due to the exponential tails of the HCS distri- 
butions these expansions are not convergent yet they are asymptotic and Borel 
resummable. To overcome this difficulty, consider the modified expansion: 

fr = fa = (^) ' e-^^'^-V^ (14) 



where 7^ = (j)a = Z^^o ^a'S'i ilau'^), and 77 > is a constant, i.e. fa is 

represented by a truncated Sonine polynomial series (whose coefficients are 

denoted by {h^^}) times a Maxwellian in which the temperature is replaced by 

a "wrong" temperature ^. It can be shown (33) that this method produces 

convergent series when 77 < 0.5. Upon substituting the form (14) in Eq. (13) 

and projecting on the Nth order Sonine polynomial S'^(7q,m^), one obtains 

2 

the following non linear algebraic system for the coefficients h'^: 



nT 



TX 



+ Xo.fi 



3mo %q 

V 3mo 



where mo = + ttlb, = and RpN is defined as follows: 



R 



'pN 



712 



dul 



■S'^i iu 



dVin 



(15) 



(16) 



10 



where 



/T^u, and the coefficients (or couphng constants), , are 



defined, in terms of the rescaled velocities u = J^^u, by: 



TjpqN 



ui2-k>0 



'a/3 



_p-vMc,ul-r)Mi3 



(iii2-k) (iui(iu2C?k (17) 



Though straightforward in principle, the computation of these coefficients 
turns out to be forbiddingly tedious as the order of the Sonine expansions 
is increased. In order to be able to carry out high order expansions, we made 
use of a computer-aided method exploiting the fact that the Sonine polyno- 
mials can be derived from their respective generating functions, Gm (■s,x): 



Gm ('5, x) 



\—m—l — -r^s 



oo 
p=0 



:i8) 



This fact enables one to define generating functions for the integrals (16) and 
(17): R{s,w) = E sPyj^RpN and B^p{w,s,t) = E s^fiw^ Elf of the co- 



p,N 



efficients, B^^f , and RpN, respectively (where all sums here and below range 
from zero to infinity). The latter can then be computed (using a symbolic pro- 
cessor such as MAPLE"^^) to any desired order by taking successive derivatives 
of R (s, w) and {w, s, t). These generating functions are given by: 



R (s, w) 



'1 — s) w 



{I - ws + [t] - I) {I - w) [l - s)f 



(19) 



and using Eq. (17) 



B. 



X 



1 , f Mo^w ^ f s \ f t \^ 



--0/3 



-lap M„ U + + 



t 



l-W 1 - S/ ' VI -t 
where lap (a, fc, c, rf, x, y, z) is the integral defined by: 

lap {a,b, c, d, x,y, z) = J (iui(iu2(ik (ui2 ■ k) 

ui2k>0 



(20) 
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where F = auj + hul + cwf + du'^ + f (ui + u'^)^ + | (ui+u'2)^ + f (ui+us)^ 
and lap (a, h) = Ia/3 {a, b, 0, 0). It is easy to calculate Ia/3, the result being: 

271-2 1 

Iaf3 = —3 7 : T (21) 

A 2 [l^afS + fJ'ap) 



where X = a -\-h + c + d -\- 2x -\- 2y + 2z, fia/s = Rap j^, i^ap = Sap 

4.aP 

— r ((^+^^ - (c + x) M^") {2Kap + t"^ ((rf+y) M"'^ - (c + x) M^^ 

i?^^ = (a + c + 2x) (Af^")^ + {b + d) [m^^)^ + i (M"^ - M^")^ (y + 2), 
Sap = - 2) (^(2c + x) (m^")^ + (M"^)^ (2c/ + y)^ + 2M"^M^°?/ 

-2 (m^")^ X , Kap = (m^° (a + c + 2x + y + ^) - (6 + d + y + 2)) , and 
= Next, using Eq. (14) for fa in Eq. (11) and (12), multiplying the 

result by s^f^ and summing over p and g, one obtains: F = y-^ "' "j_^ ^Y, 
where 



E E r^f^^B/^^B/^^, (22) 

= dsPdt^' \-'^=o,t=o, the generating function G^f^ {s,t) being given by: 

^) ^ (23) 

With the coupling constants calculated, the algebraic equations, Eq. (15), are 
now solved to produce the distribution functions, from which one can easily 
calculate the temperature ratio, Ta/Tb- Fig. 1(a) shows a plot of the tem- 
perature ratio as a function of — ^ — for a system of particles of equal mass 
density, c = 0.5 and Eaa = £bb = £ab = 0.8. Note that the temperature 
ratio can "violate" equipartition quite strongly. In order to compute the heat 
conductivity and the thermal diffusion coefficient, we now turn to solve the 
Boltzmann equation, carrying out the Chapman- Enskog expansion to first or- 
der in gradients of the hydrodynamic fields, with the homogeneous cooling 
state distribution function serving as as a zeroth order solution. The distribu- 
tion functions are expressed as = fa'^^ + fa 1 where f^ is the first order per- 
turbation in a gradient expansion, written for convenience as f^ = fa'^^ai 
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cf. Eq. (14), where the value of v may be taken to be different from that of 
T]. To first order in the gradients of the hydrodynamic fields, the two coupled 
Boltzmann equations for the perturbations read: 

where 



are the linearized Boltzmann operators, and fa denotes the first order term 
in the expansion of Dfa- Following the definition of the operator, D: 



D\nna + {r](t)a - <P'o) DVi 

3 



( 3 



DlnT 



dc 



' J a a ' a J I 



K fhLu 



D\\ina + 2'JaUil'DVi + ( ^^aU^y 



DlnT 



where a prime superscript of a function denotes a derivative with respect to its 
argument. Carrying out the expansion to first order in gradients leads, on the 
basis of tensorial considerations, to the following form for the function $^ : 



^KT f 2\ , ^^^"^ . ;^Knl \ , l9 lu n 



91nc 



'■' dxj 
dV,. 



dxj 



(24) 



where for any tensor Aij, Aij = i (^Aij + Aji - ISijAkk). D^^ is given by: 
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d\nT 



1: 

~ J a 



dxj. 



dVj I , nrrin , , , , A d\nn 



9 In c 



dc dx4 



Employing (24) in the definition of Q and J a, one obtains: 



Ka = — I r 
^ 3 



"^aJa ^a \1au ) VTA 



Uidw 



rM,u^K,T I 



wliere 3>f'"^ obey (witli P ^ a): 



dxi 



(25) 



9 JntQ ^ ^ 



{L^^l + + ^S) (7.-^) Vr.n) + (7..^) VT,u) 

/^u) = (26) 



wfiere ($f ''^^t^^u 



^2 _ nmg 



J a 



- {lau'v - I) $f v^u and 



— \4>a) u. Tlie equations (26) are solved 



subject to the constraint / /^mAU(iu+ / fj^niBudu = using the same com- 
puter aided method applied to the solution of the HCS. The functions ^^'^ 
are represented by series of Sonine polynomials ^^'^ — Ep=o ^^r—S^ [laU^) 

(the factor z/ < 1 in the function /^'^^ ensuring the convergence of the se- 
ries), which is substituted in (26). The resulting equation is then projected on 

S3, (7a'u^) -x/TaU to yield a linear system of equations for the coefficients k^'^'^. 
2 

The system of equations is then solved (after evaluating the coupling constants 
using the generating function method explained above). The resulting pertur- 
bations are then used to calculate the transport coefficients. Results for the 
thermal diffusion constant, k^, and the heat conductivity, A^, are presented 
in Fig. l(b,c). Note the nontrivial dependence of these transport coefficients 
-. Although clearly equipartition is strongly broken, these transport 



on 



mA+niB 

coefficients are well behaved (as are those not shown here). 
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Fig. 1. Plot of the temperature ratio Ta/Tb (a), thermal diffusion coefficient 
(b), and thermal conductivity A"^ (c), as a function of ^^^^^ for grains made of 
the same material, with ^ = ^ = 0.5 and eaa = ^AB = ^bb = 0.8 



4 Elastic energy and "heat flux" 



Consider a set of particles whose binary interaction potential is harmonic (or 
approximated by): U[rij{t)] = \Kij (|rjj| —Ujf', where Uj is the equilibrium 
separation of particles i and j and Yij = — r^. The force on particle i ex- 
erted by particle j is given by fj/j = —'ViU (vij) = —Kij {\rij\ — Uj) Vij, where 
= Yij/ \rij\ is a unit vector. Consider, for simplicity, an unstressed reference 
configuration, in which all particle pairs are at their equilibrium separation 
( r°- = Zjj), where the superscript denotes the reference configuration. The 
relative displacement of particles i and j is defined by: Uj^ = rjj— r?-. To linear 
order in the displacements the corresponding force can be approximated by: 



-i^,(l4 + u, 



(27) 



For atomic systems, the harmonic approximation is obtained from a Taylor 
expansion of the (effective) interatomic potential around its minimum. For 
granular materials, the contact interactions among the particles are typically 
described using contact mechanics (36; 37) (force models used in simulations 
of granular materials are reviewed in (38; 39; 40; 41; 42)). For the present 
discussion, it is assumed that the interparticle forces can be linearized for 
small deformations with respect to a reference configuration. 

Consider the contact stress, see Eq. (A. 8). Upon substituting the formula for 
the force to linear order in the displacements, Eq. (27), in Eq. (A. 8) one obtains 
an expression for the linear elastic stress tensor: 



1 



ds i 



sr 



m 
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It is interesting to examine the potential energy density as read off Eq. (A.l) 
and compare it to the relation used in continuum hnear elasticity. In the 
quasistatic hmit, the energy density reduces to the potential energy density 
[in this limit, the time of deformation goes to infinity (t/ — > oo) and the 
velocities go to zero (vj 0), so that the kinetic energy is zero while the 
displacement Uj ~ Vjtj, is kept finite]. To lowest nonvanishing order in the 

strain, the potential energy is given by \Kij (r'^j ■ u^-^ . Hence, at this order: 
e^V, ^) = ^ E (r° ■ u,,)' </.[r - r,(t)]. (29) 



Classical linear elasticity identifies the elastic energy as e^^ = : e^™, where 
e'™ is defined in Eq. (B.7). Since the stress is symmetric for central forces (it 
is easily verified that Eq. (28) is invariant to the exchange of a and /?), we 
can calculate I^^Jj^f^ = ^-£r (cq/?'^^") instead, the second equality following 

from the static equilibrium relation = 0, see also Eq. (B.3). It follows, 
using Eq. (A. 8) and definition (B.6) that: 



Next, using Eq. (B.15) one obtains: 



1 ■ ■ 1 d f 
2^SeS = E - fi/jar'jfs J ds 0[r - r° + sr° ] 



1 



Id f 

+ 4 ^ E <Ji/3o.rlp J ds 0[r - r° + sr. 



1 d 

E Uio.fil3o.r%p^^ I ds 0[r - r° + sr° ] 

1 







I „„o 1 

ij\ 



Id f 

+ 4 ^ E <afi/jarijfs j ds<P[r-v1 + sv. 



In equilibrium J2j fi/ja = 0, hence J2j fi/ja'P {'^ ~ ''^i{t)) = and it can be 
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replaced by — fi/ja'P (r ~ ^ii't))- Therefore 

= -I E ^^o.f^,a (0 [r - r,(t)] + [r - r,(t)]) 
^ ^ ij 

Id } 
+ J7^T. <Ji/jo^r^jP / ds 0[r - r° + sr°.] 

1 19 } 

= -4 E UvJi/jc.<t> [r - ri(i)] + 4^ E ^L/i/iar° /3 y rfs <p[v - r° + srj,-] 
«i ij 

Next, substituting the expression for the force, Eq. (27), one obtains: 
= i E (4 ■ u,,)' [r - r,(t)] 



+ 



4 drp 



E 



.el/ 



ds < 



sr 



E 



U/3 



(is (/)fr — 



sr,- 



(30) 



This result is identical to that obtained on the basis of dynamical considera- 
tions in Appendix B, cf. Eqs. (B.12,B.14). 

It follows that the coarse grained elastic energy density is not precisely given 
by the classical expression |cr'"^ : e^™; the above additional term provides a 
correction to the classical expression. The correction represents the adiabatic 
limit of the divergence of the (time integral of) the heat flux or, in physical 
terms, the divergence of the fluctuating part of the work of the interparticle 
forces, i.e. the work done on the fluctuating parts of the displacements or 
the work which is unresolved by macroscopic elasticity. As this term is a 
divergence of a flux, its average over a volume, Q, whose linear dimension is W, 
is proportional to W^~^, while the average over the first term is proportional to 
W^. Therefore their ratio tends to zero as 1/W, i.e., the average of |cr'™ : e'™ 
over a sufficiently large volume (not its 'local value') begets the standard 
elastic energy density. It should be emphasized that both contributions to the 
elastic energy density are of the same (second) order in the displacements (or 
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strain), i.e. the additional term is not a "higher order" correction but rather a 
correction that appears when the deformation is not affine (i.e., when 7^ 0). 



4-1 Numerical results: 2D disordered harmonic networks 

In order to examine the correction to the classical elastic energy numerically, 
we use disordered harmonic networks in 2D, obtained as follows: "particles" (of 
equal mass) are placed at the nodes of a triangular lattice, of lattice constant 
d. The system is square shaped. Nearest neighbors are connected by springs 
(with rest lengths equal to the lattice constant, ensuring a stress- free reference 
state). We introduce disorder in two ways: first, the spring network is diluted 
by randomly removing a fraction c of the springs. Second, we employ spring 
constants which are uniformly distributed in the range [K — SK, K + SK]. The 
particle positions are the same as in the triangular lattice, hence the density 
remains uniform. We use periodic boundary conditions in both axes, and apply 
a specified global strain to the system by imposing a change in the period as 
well as (for shear deformation) Lees-Edwards boundary conditions. 

The (linear) static equilibrium equations are solved (by inversion of the dy- 
namical matrix) for a given applied global strain, e'^fj, yielding a set of cor- 
responding displacements {uj}. The latter are used for calculating the CG 
linear strain field, e fusing Eq. (B.7)], the CG linear stress, a^'"" [Eq. (28], 
the CG potential energy density, e^' [Eq. (29)], and the fluctuating energy 
flux, iQfoi^ccrcv,! [j^q^ (B.14)]. The coarse graining function used (43) is a 2D 

Gaussian, (/)(r) = ^e-^'"!/'^)'. 

We verified that for a lattice with uniform spring constants (i.e., c = 6K = 0), 
the strain components are uniform and equal to the applied global strain for 
w > d, a.s expected, since in for this system, under homogeneous deformation, 
the particle displacements are affine. The resulting stress components are uni- 
form as well, and are consistent with the continuum isotropic elastic moduli for 
a triangular lattice (Lame constants X = fx = -^^-|^). As expected, in this case 
the CG energy density conforms to the classical expression, = |cr'"^ : e^™. 

For disordered systems, even under a "homogeneous applied strain" as de- 
scribed above, there is a correction to the classical expression for the en- 
ergy density, as detailed above. Fig. 2(a) presents a contour plot of the rel- 
ative energy difference Ae = 2 ^Hn^gun — 1 = — '^"^un.^iin ' ; for a system of 
30 X 30 particles, with c = 2% and 6K/K = 25% subject to an applied strain 
e^PP = -e^PP = 5 ■ 10-3, e^PP = 0, calculated with a CG width, w = Ad. The de- 
viation from the classical expression is non- negligible (Ae ~ 1%), even thought 
the present system has a relatively small disorder [we verified that in the same 
system, a linear stress-strain relation holds to within 4% (43)]. It is interesting 
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Fig. 2. (a) A contour plot of the relative difference between the elastic energy density, 
e*^', and its "classical" value, ^c^™ : e^™, in a disordered harmonic network (see text), 
(b) The fluctuating energy flux, lQforce,rev,i (bi^ck) and the (negative) gradient of 
the "noise", — Vr/ (gray), in the same system (arbitrary scaling). 

to try to relate the fluctuating energy flux to a local measure of the disorder in 
the particle displacements. A possible field for characterizing the fluctuating 
displacements is the "noise" field (44), ?7(r) = X)j '^j|u'j(r)p0[r — rj(t)], defined 
in an analogous way to the kinetic temperature. In Fig. 2(b) we compare the 
fluctuating energy flux with the gradient of the noise field. It appears that 
their magnitude (but not their directions) exhibits quite a strong correlation, 
suggesting that the fluctuating energy flux may indeed be related to gradients 
of the noise, but the corresponding "conductivity" may be anisotropic. This 
is to be expected, since the system is locally anisotropic [as the local elastic 
tensor has been checked to be anisotropic (43)]. It is possible that a tenso- 
rial, rather than scalar, characterization of the displacement fluctuations is 
required, analogous to an "anisotropic temperature". 



5 Conclusion 

We have shown that the classical definition of granular temperature as a mea- 
sure of the fluctuating kinetic energy is useful for the description of granular 
gases even when this temperature is not the (inverse of the) thermodynamic 
conjugate potential of the energy, as in (local) equilibrium. An analysis of the 
case of a near-elastic granular gas shows that indeed one can venture into the 
nonequilibrium domain starting from local equilibrium (for the near elastic 
case). The example of a binary granular gas mixture in which the tempera- 
tures of the components can be very different from each other and still the 
Chapman- Enskog method produces sensible results, is a compelling case in 
favor of the use of the granular temperature field for far-from-equilibrium sys- 
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terns, at least in the fluidized phase. In the 'opposite' hmit of a granular solid 
we have shown that the heat flux (which is always present, even in athermal 
system, as shown in Appendix A) is not 'lost' when one performs a defor- 
mation of an elastic solid: its effect appears in the form of a correction to 
the classical elastic energy density. This correction is strongly correlated with 
the gradient of an 'elastic temperature' which is a measure of the degree of 
non-affinity induced by the randomness of the system. Note that the system 
considered is not pre-stressed whereas many granular systems are. This fact 
does not change the fundamental conclusion that in elastic solids there is an 
accumulation of energy as a result of work done on the non-affine (or fluctu- 
ating) degrees of freedom. As the physics of this flnding is very general we 
believe it applies (with a proper modiflcation) to granular solids beyond the 
elastic regime. Indeed a strong correlation of the displacement 'noise' with the 
appearance of localized plastic events has been suggested in (44). Another con- 
nection of this kind has been proposed in (45). At this stage we do not know 
whether there is a direct link between this 'new' contribution that stems from 
the heat flux and (gradients of) 'conflgurational temperatures'. However, since 
there is an obvious link between diffusion and non-affinity (see e.g., (46)) and 
a seemingly established link between diffusion and a 'conflgurational tempera- 
ture' (through an FD relation) it comes to reason that the 'extra' term in the 
elastic energy is related to a conflgurational measure of disorder. Following the 
spirit of a suggestion by Savage (47) one may perhaps conjecture that an effec- 
tive temperature which combines the conflgurational and kinetic contributions 
(perhaps just as a sum) can serve as a generalized granular temperature valid 
for all granular phases (with appropriate extensions to other systems). At this 
stage it seems that much more research is necessary to establish whether this 
or similar conjectures are of any physical value. As a flnal comment we would 
like to mention the rather well known fact that (even relatively simple) non- 
equilibrium systems exhibit rich behavior and cannot usually be characterized 
by the same number of macroscopic variables or flelds as (near) equilibrium 
systems. Therefore, even if it turns out that many properties of granular mat- 
ter can be correlated with a generalized temperature or temperatures there 
is a good chance that more characterizations will be needed for this class of 
systems to fully describe or explain their behavior. 



A Continuum Mechanics: a microscopic derivation 

The equations of continuum mechanics are derived below from microscopic 
considerations. For sake of simplicity it is assumed that the constituents expe- 
rience binary interactions whose reversible part is determined by a potential. 
The derivation below is a slight extension of what can be found in (26; 48) 
but unlike in (26) we do not invoke temporal coarse graining, for simplicity. 
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Consider a set of particles of mass rrii, whose (center of mass) positions and 
velocities are rj(t) and Vj(t) = rj(t), respectively. Let the binary interaction 
potential be denoted by U{rij), where Tij = rj — r^. It is assumed that U = 
for i = j. In addition, let 0(r) denote a spatial coarse graining function, i.e. is 
normalized (its integral over space equals unity), it is positive semidefinite and 
has a single maximum at r = 0. It is often convenient to employ a Gaussian 
distribution for </>. The width of (p is the "coarse graining scale". 

The coarse-grained (CG) mass density, p(r, t), at position r and time t is given 
by: p(r, t) = J2i ^i'P [r ~ ^i{t)\- Note that while rj(t) is the position of particle 
i at time t, r is merely a considered position in space. Similarly, define the 
CG momentum density, p(r, t), by: p(r, t) = I]j mjVj(t)0 [r — rj(t)]. The CG 
energy density, e(r,t) is given by: 

e(r, t)^Y. \^^^^<P [r - r.(t)] + ^ E U{v,,{m [r - r.(t)] (A.l) 

By taking the time derivative of the density one obtains: 
d 

p{r,t) = -—Y,miVia{t)<l)[r - r,{t)] = -divp(r,t) (A.2) 



where the summation convention has been invoked and Greek letters denote 
Cartesian coordinates. Defining the velocity field, V(r, t), as: V(r, t) = , 
one obtains from Eq. (A.2) the equation of continuity: p = — div(pV). Next, 
upon taking the derivative of the momentum field one obtains: 

d 

pair, t) = -— E miVia{t)vip{t)(j) [r - Yi{t)] + ^ miVi^cj) [r - Yi{t)] (A.3) 

i i 



Denote by fj(t) the force experienced by particle i and by ii/jit) the force 
exerted by particle j on i. Using Newton's second and third laws one obtains: 
Ei miVi^{t)(f) [r - Tiit)] = fi/ja{t)(j) [r - r.i{t)] = fj/ia{t)(p [r - rj(t)] = 
- Eij fi/ja{t)4> [r - ri(i)]- It follows that: 

^miVia{t)(t) [r - ri{t)] = 



= 7;i: h/Ut) (0 [r - r,(t)] - [r - r,(t)]) 

1 } d 

= - o E fi/joit) / ds—J) [r - nit) + SYijit)] 
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d 1 



± 

Yl ^ifi/jarijpit) J ds(f) [r - ri{t) + srij{t)] (A.4) 



Substituting the latter result in Eq. (A. 3) one obtains: Pa(r, t) = — ^Pa/3(i", t), 
where the pressure tensor, P, is given by: 

P„/3(r,t) = 'YmiViait)vip(t)(t> [r - ri(t)] 



1 } 
+ 2 ? ^) J ds(f) [r - ri{t) + srij(t)] (A.5) 

ij 

Next, define the fiuctuating part of the velocity of a particle (when coarse 
graining is around the position r at time t) as: v-(r,t) = Vj(t) — V(r, t). 
Substituting this definition in Eq. (A.5) and performing straightforward alge- 
bra (using I]i ^^(r, t)0 (r — rj(t)) = 0, following the definitions of v' and the 
velocity field V), one obtains: 

Pap{r, t) = p(r, t)\4(r, t)Vp{r, t) + ^i^Li^, t)v'if,{r, t)0 [r - r,(t)] 

i 

1 } 

+ 2 ^ tyijfsir, t) J ds(f) [r - ri{t) + srij{t)] (A.6) 

ij 

Identifying 

^ap{r,t) = - Y.mvL{r,t)v:^{r,t)<P[r -r^it)] (A.7) 

i 

as the "kinetic stress" and 

1 } 
aX(r, ^) = -2 E t^i^r, t) J ds<P [r - r,(t) + sr,,(t)] (A.8) 



as the "contact stress", the stress tensor, cr, is given by: 
a„/3(r,t) = -^mii;,'„(r,t)t;,'^(r,t)0[r-ri(t)] 

i 

1 } 
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Note that an opposite convention for the sign of the stress tensor is often 
employed in the field of granular matter. It follows from Eq. (A. 9) (employing 
p = pV and invoking the equation of continuity) that: 



P 



Dt 



(A. 10) 



where = ^ + V ■ VV denotes the material derivative. 

Next, consider the energy density. Taking the derivative of Eq. (A.l) one 
obtains after some simple algebra and the use of Newton's second law: 



e(r,t) 



d 



dr. 



E 



Vic,{t)(j)[Y - Yi{t)] 



Y: m ■ Mm [r - r,{t)] - - E il/Jit) ■ v.,(t)0 [r - r,(t)] (A.ll 



where the reversible part of the force fi/j is given by fjyj(t) = — ^ • Using 

the above decomposition of the velocity, Vj(t) = v^(r, t) + V(r, t), one obtains 
that the term whose divergence appears in Eq. (A.ll) can be rewritten as: 



^ ^m,v^{t)vUt)^ [r - r,(t)] + U{r,,{t))v,^{t)4> [r - r,(t)] 

i ij 

= e(r, t) V,(r, t) + ^ E ^^^'(r, ^ [r - r,(t)] 



- Vfsir, t)a^-(r, t) + - E U{r,,{t))vUr, t)<j> [r - r,(t)] 
Next, note that 

y: m ■ Mm - ^^it)] = e %(^) ■ Mm - Mt)] 

i ij 

= -E%(t)-v,(t)0[r-r,(t)] 



(A.12) 



5:f./,(t)-v,,(t)0[r-r,(t)] 



+ 7 E%(^) ■ (Mt) + v,(t)) (0 [r - r,(t)] - [r - r,(t)]) 
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^^%(t)-v,,(t)0[r-r,(t)] 



d 1 
9r„ 4 



± 

^ ■ [v,(t) + v,-(t)] r,,,(t) 1 [r - ri(t) + sr^,{t)] (A. 13) 

«i 



Using velocity decomposition and Eq. (A. 8) one obtains from Eq. (A. 13): 
Y: m ■ v^(t)0 [r - r,(t)] = W %(t) ■ v,,(t)0 [r - r,-(t)] 

i ij 



1 d 

4 9rQ, 



d 



%(^) ■ [v-(r, t) + v;.(r, t)] ri^„(t) j ds(f) [r - ri(t) + sv,^{t)] 



+ ^(V/3(r,t)a--(r,t) 



(A. 14) 



Next, combining Eqs. (A. 9, A. 11, A. 14) one obtains: 



e(r,t) 



d 



e(r, t)V,(r, t) - V^(r, t)a/j„(r, t) + g„(r, t) - p(r, t)r(r, t)(A.15) 



where Q is identified as the heat flux and is given by: 

+ 7E%W- (v^(r,^) + v;(r,t))r,,„(t) j ds0 [r - r,(t) + sr,,(t)] (A.16) 

ij 

and the sink term, —pT is given by: 

-p(r,t)r(r,t) = (%(r,t) -f-^(r,t)) ■ v,,(r,t)</.[r-r,(t)] 

1 



E^-^(r,t)-v,,(r,t)0[r-r,(t)] 



(A. 17) 



where the dissipative part of the force, f™^, is given by fj/j — f[j^^. Note that 
the heat flux represents the transport of energy by the fluctuating velocity 
and the power related to the fluctuating velocities. The energy sink term is 
entirely due to the power associated with the dissipative part of the forces. 
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The energy density can be rewritten as follows: 

i 

+ ^p(r)V^(r, t) + i ^ f/(r.,(t))0 [r - r.(t)] (A.18) 



Denote the internal specific energy by u. The latter is given by 

p(r, t)n(r, t) - + ^ E ^^(r., [r - r.(t)] (A.19) 



and the energy density therefore equals |pV^ + pw. Substituting this decom- 
position into Eq. (A. 15) and using Eqs. (A. 2, A. 10) one obtains: 

P(r,t)^^^^^ = ^^|^a/3«(r,t) -divQ(r,t) -p(r,t)r(r,t) (A.20) 



B The quasistatic limit of Continuum Mechanics 



Consider the case of slow deformations of a near-static assembly of particles. 
To this end, let the particle velocities be denoted by Vj(t) = 5vj(t), where 5 is 
a small parameter and A is considered to be 0(1) for any dynamical variable, 
A. The velocity field, V, as well as all other fields are rescaled in a similar 
way. Furthermore, let time be expressed as t = |, hence ^ = 5-^. With this 
rescaling the equations of continuum mechanics, Eqs. (A.2,A.10,A.15) become: 

5^ = -Miv (pV) 

^2 DV^ ^,dal'^ , dalj' 
dp — ^ ' 



Dt drp drp 



dt dt dfc 

S^fi'^pT + SQl"' + + SQ'r'' ) - 5pT (B.l) 



where the superscript 'kin' refers to the kinetic part, the superscript 'pot' 
refers to the potential part (involving the potential U) and the superscript 
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force refers to the fluctuating power contribution to the heat flux. To linear 
order in 6 the above equations reduce to: 

^ = -div (pV) (B.2) 

n cont 

"^^-^ -0 (B.3) 



^gPOt Q 



,poty^ - V^ga^r + QT + Q'rn - Pf (B.4) 



dt dvc 



Using Eq. (B.3), the equation for the energy density becomes to linear order 
in the 'slowness', 5 (and after reverting to the original time variable): 



n„pot rl / \ ri\f 



We further specialize to the case when the particle relative displacements are 
small and the potential is quadratic in the relative displacements (to leading 
order). Let Uj(t) be the displacement of particle i: Uj(t) = Vj. Also, define the 
displacement field u(r, t) via the Lagrangian relation: ^"'^f'^^ = VL-(R,t), 
where R is the initial position of a pathline and the superscript 'La' denotes 
the fact that one is employing Lagrangian coordinates. As shown in (43) 

„nn( )^ S.m.u.(t)0[r-r.(t)] 



where the superscript 'lin' denotes the linear order in the displacements. To 
this order the (linear) strain field, e'"^ is given by: 



Keeping in Eq. (B.5) only terms that are quadratic in the displacements one 
obtains (note that the potential is of second order in the displacements): 

= -div Qf- + ^a™- - pr (B.8) 

Since for quadratic potentials cr^°^^ is symmetric (as is easy to check, see also 
Eq. (28)) one can rewrite Eq. (B.8) as follows: 



dt 2 V dra drp 
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hence (to quadratic order in the displacements) 

n pot ^^lin 

^ = -div + -^af:' - Pr (B.IO) 



Consider next the following virtual dynamical evolution from a force free initial 
state to a final state at time tf. Uj(t) = Uj(tf)^. Assume that the reversible part 
of the force, f^yj is linear in the displacement difference Uj^ and the irreversible 
part of the force is linear in the velocity difference Vjj. In this dynamics the 
velocities are constant, the strain field satisfies e'™(r,t) = e^™(r, tf)^, the 
reversible part of the contact stress, cr™"*'''^^ (in which the force is replaced by 
the reversible part of the force) satisfies cr'^°°*'''''^(r, t) = cr™'^*'''°^(r, tf)^, the 
irreversible part of the contact stress, cr™'^*'"'''°^(r, t), is constant in time, and 
similarly the part of the heat flux which corresponds to the reversible parts 
of the forces satisfies Qf°'"^'=''''=^(r, t) = Q^°'"'="''^°^(r, tf)^ and Qf°'"=^'™^(r, t) is 
constant in time. In addition, pr(r, t) is constant in time. Note that within the 
approximation in which only the lowest nonvanishing order in the displacement 
is retained, rjj(t) is replaced by rjj(O). Upon substituting these expressions in 
Eq. (B.IO) and integrating Eq. (B.IO) over time from zero to tf one obtains: 



/ I \ lin / 4. \ cont,rev/ ±\, lin / ±\ cont,irrcv/ • \ 

e(r, U) = 2^a/3(''' ^f)^/3a (r, ^f) + ea/3(r, ti)(T^^ (r, ti) 



^div Qfo-.-(r, tf)tf - div Qfo-.--(r, t,)ti - p(r, t,)T{v, t,)t, (B.ll) 



In order to keep the displacements fixed we replace the velocities that ap- 
pear in the expressions for the various fields by l/U times the correspond- 
ing displacements. A field in which the velocities have been replaced by the 
corresponding displacements is marked by a superscript I. It follows that at 

j. ^cont.rcvjl ^cont,rcv ^cont,irrGv,I ^^^cont.rev Qforce,rev,I ^^QforcG,rcv 

Qforce,irrev,i ^ ^2Qforce,irrev (j.Qf.r^\\ the assumptlou that fi^J^^^ is linear in the par- 
ticle velocities), and (pE)^ = tf pE. Substituting these expressions in Eq. (B.ll), 
keeping the displacements fixed and letting tf go to infinity, one obtains: 

e(r, tf) = U^,{v, U)alT"'\^. U) - ^div Qf— .-.^(r, t,) (B.12) 



where (changing tf to t) if follows from the expression for the heat flux, 
Eq (A. 16) and the above definitions that 
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1 

X J ds<t) [r - r,(t) + svij{t)] (B.14) 



where 

u:(r,t)^u,(t)-u'-(r,t) (B.15) 



is the displacement fluctuation of particle i. We have thus obtained a correction 
to the classical formula for the elastic energy which represents the fluctuating 
work and stems from the integral of the heat flux over time. This result can 
be directly obtained from an analysis of the elastic energy, cf. Eq. (30). 
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